function [v_ratio,tot_var] = var_decomp(irfa)
% PURPOSE: Variance decompositions
% Computes the Forecast Error Variance Decomposition (FEVD) 
% for orthogonal shocks.
% -----------------------------------------------------
% USAGE: [v_ratio] = var_decomp(irfa)
% where:
% IRFA = Input is the Impulse Response Function for the system.  This can
% be obtained from the MK_IRF procedure.
% -----------------------------------------------------
% RETURNS: 
% This proc returns "v_ratio", the FEVD as an (K*K)x(MALAGS) matrix. This is MALAGS blocks.
% The i-th row of block h-th is associated with the h-step ahead forecast of y_i, 
% the j-th column is associated with the j-th impulse.
% sum(v_ratio,2)

H=size(irfa,3); %Horizon
kk=size(irfa,2);

%Column j contains the squared IRFA of variable j
irfa_sq = reshape(irfa.^2,kk,kk*H);

v_ratio=zeros(kk,kk,H);

d=zeros(kk,kk);
tot_var=zeros(kk,H);
 
for t=1:H
    %d CONTAINS THE PARTIAL VARIANCE DUE TO SHOCK 1 (column 1), SHOCK 2 (column 2), and so on
    %of variables 1 (row 1), of variables 2 (row 2), and so on.  
    %[(THETA11,1 + THETA11,2 +...+ THETA11,t), (THETA12,1 + THETA12,2 +...+ THETA12,t);
    % (THETA21,1 + THETA21,2 +...+ THETA21,t), (THETA22,1 + THETA22,2 +...+ THETA22,t) ]
    d=d+irfa_sq(:,(t-1)*kk+1:t*kk);
    
    %Matrix TOT_VAR contains the total variance up to point 'h' in column 'h'
    %of variables 1 (row 1), of variables 2 (row 2), and so on.  
    %(THETA11,1 + THETA11,2 +...+ THETA11,h + THETA12,1 + THETA12,2 +...+ THETA12,h);
    %(THETA21,1 + THETA21,2 +...+ THETA21,h + THETA22,1 + THETA22,2 +...+ THETA22,h)
    tot_var(:,t)=sum(d,2);
    
    for ii=1:kk
        v_ratio(ii,:,t)=100*d(ii,:)./tot_var(ii,t);
    end
end
